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Abstract 

We study gravitational collapse of a charged black hole in f(R) gravity using double-null 
formalism. We require cosmological stability to f(R) models; we used the Starobinsky model 
and the R + (l/2)cR 2 model. Charged black holes in f(R) gravity can have a new type of 
singularity due to higher curvature corrections, the so-called /(i?)-induced singularity, although 
it is highly model-dependent. As the advanced time increases, the internal structure will 
approach the Cauchy horizon, which may not be an inner apparent horizon. 

There is mass inflation as one approaches the Cauchy horizon and hence the Cauchy horizon 
may be a curvature singularity with nonzero area. However, the Ricci scalar is finite for an out- 
going null observer. This can be integrated as follows: Cosmologically stable higher curvature 
corrections of the Ricci scalar made it bounded even in the presence of mass inflation. 

Finally, we conjecture that if there is a general action including general higher curvature 
corrections with cosmological stability, then the corrections can make all curvature components 
finite even in the presence of mass inflation. This might help us to resolve the problem of inner 
horizon instability of regular black hole models. 
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1 Introduction 



Quantization of gravity is the final goal of modern physics. To quantize gravity, there exists two 
non-trivial difficult problems: gravity is not renormalizable and gravity allows various singularities. 
In this paper, we try to combine these two problems in an approximated way. 

Traditionally, it is well-known that gravity is difficult to quantize due to the problem of renor- 
malization. When we calculate scattering amplitudes for a given action, we have to draw and 
calculate Feynman diagrams. If a Feynman diagram has a 'loop', then the contribution of the 
diagram simply diverges. By renormalizing parameters suitably, we can cancel infinities and can 
fit with the observable scattering amplitudes. The problem of quantization of gravity lies in the 
very fact that, gravity is non-renormalizable. This implies that if the related energy scale is greater 
than the Planck scale, and hence if the effect of quantum physics and gravity are of similar orders 
of magnitude, the effective action then has, not only finite terms, but also infinitely many terms to 
cancel infinities up to leading order. Therefore, such a theory cannot give observable results when 
the energy scale becomes greater than the Planck scale. 

The problem of renormalization may be resolved if we accept string theory. It is an well- 
established fact that scattering amplitudes of string theory is naturally finite up to all orders. In 
the low energy limit, string theory gives effective action in the form of supergravity. If we want to 
know higher order terms of the low energy effective action, we can perturbatively expand it and 
evaluate all coupling parameters from string theoretical calculations in principle. For example, the 
effective action of the gravity sector will have the following form [1] : 

S = 5Ap/^V^-*[(*+(V0 a ) 
+ a' (aiRl„ p<T + a 2 R% + a 3 R 2 

+ ai R^V p <t>V v <t> + a b R{V4>) 2 + a 6 i?V 2 + a 7 (V 2 0) 2 + a 8 V 2 0(V0) 2 + a 9 (V(/>) 4 ) 
+ a' 2 0(R 3 )] 

+ [R+...]+e+[R + ...] + ... , (I) 

where A s is the string length, a' is the string tension parameter, d is the space dimension, <f> is 
the dilaton, and a;S are numerical coefficients. Here, what we want to say is that it is too difficult 
to deal with this general problem without any approximation. Therefore, we have to reduce this 
action for meaningful calculations. 

Let us assume d = 3 for simplicity. Any reasonable model of string theory should be reduced to a 
four-dimensional model. In this paper, we will not consider dilaton and other sectors for simplicity; 
we will only consider the gravity sector including a U(l) gauge field and a complex scalar field 
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to study a charged black hole. We have to choose how to include the higher curvature terms. 
There can be infinitely many combinations of higher curvature terms: R, A = R^ pal B = Rz u , 
C = RftvpaR^ 1 ' R pa , etc. The action should have the following general form: 

S= ^kJ dx '^[I{R) + h(A) + f 3 (B) + U(C) + ...], (2) 

where fiS are arbitrary functions. 

At this point, it is important to stress on some general results of the f(R) gravity. If we choose 
only the f(R) term, then the above action can be changed to that of the Brans-Dicke [2] type by 
a simple field redefinition (3J. Here, the Brans-Dicke field is related to a non-trivial potential and 
this potential has further cosmological implications. To make the model cosmologically viable, we 
have to impose some conditions on the function / [3]. Imposing such constraints, it is well-known 
that f(R) gravity can be a viable model of cosmology up to the status of current observations. 

Of course, there exist many difficult theoretical problems in this context. Some of them are as 
follows: (1) What will be the back reactions of higher curvature terms, which were neglected by 
f(R) gravity? (2) Is the model really renormalizable? (3) Can the model be embedded into string 
theory? They are difficult to answer. However, at least we can say the two things for general viable 
f(R) models as follows: 

• Introducing an auxiliary field, we can redefine the gravity sector to be (non-minimally coupled 
scalar field) + (Ricci scalar) + (non-trivial potential of the auxiliary field). 

• The non-trivial potential should allow a stable minimum to satisfy cosmological constraints. 

It is not known whether we can generalize these two things for the general action (Equation ©), 
although we think that it will not be too exotic generalization. 

In this paper, we argue that the two conditions for f(R) gravity is useful to resolve a type of 
singularity - the mass inflation singularity. From the singularity theorem, we know that, if we 
assume a reasonable causality condition (e.g., global hyperbolicity), a reasonable energy condition 
(e.g., the null energy condition), and a reasonable initial condition (e.g., a trapping horizon), then 
there should be a singularity [5]. There are some known ways to resolve such singularities; violating 
some of the assumptions of the singularity theorem, or assuming a certain quantum gravitational 
correction. These models are known as regular black holes 6 7 . However, one of the difficult 
problem in case of regular black hole models is, not only do they have the outer apparent horizon, 
but also, they do have an inner apparent horizon. The introduction of the inner apparent horizon 
is useful in resolving the problem of central singularity. However, the inner horizon is in general 
unstable via mass inflation [8] [9]. In other words, some of the curvature functions diverge at the 
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inner horizon, and then the model easily goes beyond the range of validity (e.g., Planck scale) of 
the initial theory. As it stands, there is no known mechanism to resolve such a mass inflation 
singularity. Perhaps, in case of large number of massless fields, we can trust this model beyond 
the inner horizon [10] . However, for a realistic application, it is difficult. This might imply that, a 
regular black hole 'metric' description always fails, the central region of a black hole is fuzzy and 
hence we need full quantum gravitational description [11) . 

However, higher curvature corrections may shed some light to the problem. Higher curvature 
corrections may be helpful in regularizing the mass inflation singularity. In this paper, we assume 
cosmologically stable f(R) gravity models and study gravitational collapse of a charged black hole. 
There is mass inflation and some of the curvature components may diverge. However, we observe 
that the Ricci scalar becomes finite even near the Cauchy horizon. This feature is derived from the 
two properties: we include higher order curvature corrections of the Ricci scalar and we impose 
the cosmological stability. Perhaps, we may generalize that general higher curvature corrections 
to a cosmologically stable universe can make all orders of curvatures to be finite near the Cauchy 
horizon of a charged black hole. Then, the initial action can be still valid even in the presence of 
mass inflation. This conjecture will be helpful to understand the mass inflation singularity. 

This paper is organized as follows: In Section [51 we describe a model for f(R) gravity and im- 
plement the model to double- null formalism to study gravitational collapse. In Section [31 we report 
gravitational collapse in f(R) gravity and discuss mass inflation, and in Section [4l we interpret and 
discuss our results. 



2 Double-null formalism for f(R) gravity 



2.1 f(R) gravity 



The action of f(R) gravity can be written as 




(3) 



Introducing an auxiliary field f/ 1 , we can change the gravity sector by 




(4) 



with a constraint ip = R. If we define a new field $ by 



* = /'WO, 



(5) 



then we can rewrite the action as 




(6) 
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where 



y(<fr) = -/(v>) + V>/'(# (7) 

This is exactly the ui = limit of the Brans-Dicke theory. 
2.2 Brans-Dicke theory 

The action of the Brans-Dicke theory with a scalar field becomes 

Sbd = J dx^J—g ($R - |*„ i * ; „«r - V(*j) 



(8) 



where R is the Ricci scalar, <E> is the Brans-Dicke field, and is a complex scalar field with a gauge 
coupling e and a gauge field A^, where F^ v = A v -^ — A^ v . Here, lo is the Brans-Dicke coupling 
constant and we choose lo = 0. 

The Einstein equation is as follows: 

G>l/ - SnT^ + 8tt-^ = 8nT^, (9) 

where the Brans-Dicke part of the energy-momentum tensor is 



iBD 



(10) 



a»" 87r$ y^';pff.y 1 1 87r< j,2 I 2 - yAi ;p ; - y 7 yp 167r 

and the matter part of the energy-momentum tensor is 

+ \ {-4>;nizA v (j) + (fr.JeA^ + 4>.^ieA„(j) - ^ieA^tp) 
+ -t F »P F » P + e 2 A^# + £™9^- (11) 
The field equations are as follows: 

^-J C -J^(^- 2V ) ~ °' (12) 

(j>;„ v g^ + ieA^(2^ + ieA^)+ieA^g^4> = 0, (13) 
7^F V KV - iecj) - ieA^J) + iej + ieA„0) = 0, (14) 



where 



T c = T c ^. (15) 
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2.3 Choice of potential V($) 

In order to make a f(R) model free from ghost, we consider $ = f'(R) > 0. Additionally, for 
the consistency of the field redefinition, we only consider a region f"(R) > 0. Dynamics of the 
Brans-Dicke field is governed by $V'($) - 2V($) = 2f(jp) - ipf'(ip). In other words, if T c = 
and $V($) — 2F($) = 2f(iji) — ipf'(ip) = 0, there can be an extreme point or a constant curvature 
point R = Rq. For convenience, we define an effective potential U($>) as 

£/($)=/ ($V" (*) - 2V(*)) d$. (16) 

Then the equation of the Brans-Dicke field becomes: 

$ ^ _ _^ T C - -J_^ = (17) 



and the stability around the extreme point becomes 

d 2 U ff'W 



d<P 2 



1> 



> 0, (18) 

ipo 



where <J> and ip are solutions of $V($) - 2V($) = 2/(^) - 0/'(V>) = 0. 

2.4 Implementation of the double-null formalism 

We use the double-null coordinates 

ds 2 = -a 2 (u, v)dudv + r 2 (u, v)dfl 2 , (19) 

assuming spherical symmetry. Here, u is the retarded time, v is the advanced time, and 9 and ip 
are the angular coordinates. 

We follow the notation of previous papers [H] [E] [E] [H] [S] : the metric function a, the radial 
function r, the Brans-Dicke field $, and a scalar field s = y/incf), and define 

Q/ Q/ 

h=—, d=—, f = r. u , g = r tV , W = $, u , Z = $ tV , w = s^ u , z = s v . (20) 

The Einstein tensors are then expressed as follows: 

Gun = -~{f,u-2fh), (21) 
r 

Guv = ^(4rf, v +a 2 +4fg), (22) 

G vv = --(g. v -2gd), (23) 
r 

Gee = -4^(d, u + Zz). (24) 
or \ r I 
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Also, we can obtain the energy-momentum tensor for the Brans-Dicke field part and the scalar 
field part (co = 0): 



^BD 



T, 



BD 



87T$ 



(W, u - 2hW), 



Z, u gW + .fZ + a 2 V 



47rr$ 



8tt$ 

{Z tV -2dZ), 



32tt$' 



8tt$ 

J2 



27ra 2 $ 



47ra 2 $ 



(</W + /Z) 



16tT$' 



(25) 
(26) 
(27) 
(28) 



rp<J 



rpC 



\ww + iea(ws — ws) + e 2 a 2 ss\ , 



47T 

Ana 2 ' 



Ana 2 



2f a j 

(hi! + zw) + iea(zs — zs) H f=?— 



(29) 
(30) 
(31) 
(32) 



To implement the double-null formalism into the numerical scheme, it is convenient to represent 
all the equations as first order differential equations. Note that 



1 ~ Q ,2 i ™ + r 2 J e» 



(33) 



and 



(34) 



The Einstein equations for a. 



and the field equation for $ are then coupled: 



/ 1 


l/r 


1/$ \ 


' da 


\ 




(* 


\ 





1 


r/2$ 












V 





r ) 


V Z -u 


) 






J 



(35) 



where 



21 



27ra 2 r 1 l, m a 2 V 



47rr, 



-T 



/5 1 



a 

4r r 
47rr / a 2 



^W + fZ) + —V, 



-fZ-9W-^(^TC + — (ZV>-2V) 



(36) 
(37) 
(38) 
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After solving these coupled system of equations, we can write all the equations as follows: 

(39) 





f.u 


= 2fh 




g,v 


= 2gd 




= h v 


= a- 


g,u 


= f,v 


= 93- 




= w lV 


r 



-^(W, u -2hW)-—T- u 

23 

€_ 

2$' 



including the scalar field equation 



a 2 q 

„2 



(40) 

(41) 
(42) 
(43) 

(44) 



Q,v = --^-{sz-sz), (45) 

fz gw iearz ieags ie 2 

z, u = to,* = 1 ■ T^a gs- (46) 

f 7* f 7* Z^y^ 

In fact, for a given f(R), we have to find the inverse function of f'(R) to write as a function 

of $. In general, it is not convenient for a numerical calculation. Therefore, we will change terms 
depending on V($) with the help of the following identities: 

V($) = -f(R) + Rf'(R), (47) 

^-2V = 2f(R)-Rf'(R), (48) 



where the Ricci scalar in the double-null coordinates can be expressed as 

V a J & r 



R = X 



(49) 



In addition to this, we use the following two equations to control the evolution of Ricci scalar R: 

= 7wy (50) 

R, = $ty (3D 

Now, the equations for a iUV , r iUV , $ iU „, and s jUV parts can be represented as first order differential 
equations. We can then implement the same integration scheme as was used in previous papers 
[12] [15] [I3| |I5] to solve the Brans-Dicke theory. We use the second order Runge-Kutta method to 
solve the problem |17j . Tests of the convergence are provided in the Appendix. 

2.5 Initial conditions and free parameters 

We need initial conditions for all the functions (a, h, d, r, /, <?, $, W, Z, s, w, z, a, q, R) at the initial 
u = U\ and v = v\ surfaces, where we set U{ = V\ = 0. 
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We have the gauge freedom to choose the initial r function. Although all the constant u and 
v lines are null, the freedom to choose the distances between them still remains. Here, we choose 
r(0, 0) = ro, f(u, 0) = r u o, and g(0,v) — r v o, where r U Q < and r„o > such that the radial 
function for an in-going observer decreases and that for an out-going observer increases. 

First, we assume that the Brans-Dicke field is asymptotically $o, where $qV'(§o) — 2V($o) = 
2f(R ) - R f(R ) = 0. Then, $(u,0) = $(0,«) = $ , W(u,0) = Z(Q,v) = 0, and R(u,0) = 
i?(0, v) = Ro- We use a shell-shaped scalar field, and hence its interior is not affected by the shell. 
Thus, we can simply choose s(u, 0) = 0. Also, w(u, 0) = h(u, 0) = a(u, 0) = q(a, 0) = holds. Since 
the asymptotic mass function 

m(u,v) s l(l + 4^ + ^-^r») (52) 

should vanish at u — v — 0, it is convenient to choose r„o — —1/2 and r v o — 1/2, and then 

We need more information to determine d, g, Z, and z at the v = surface. We obtain d from 
Equation pTj) . g from Equation (fl2"j) . Z from Equation (|43[> . and z from Equation (jlS)) . 

We can choose an arbitrary function for s(0,i>) to induce a collapsing pulse. In this paper, we 

use 

fo v ~ Vi \ , ■ ( i v ~ Wi 

COS 27T + Z COS 27T- 



s(ui, v) = A sin 7r 



Uf - Vi 



(54) 

Vf — Vi J \ v f — v i 

for < v < V{ and s(0, u) = otherwise, where V{ is the width of the pulse, A is the amplitude, 
and 5 is a free parameter to tune the phase of the complex scalar field [13] [II]- Then we obtain 
z{ui,v). This implements one pulse of energy (T vv ~ z 2 ) along the out-going null direction by a 
differentiable function z(0,v). Also, from Equation pi)]) , we can use d = r\z\ 2 /2g& on the u = 
surface and thus we obtain d(0,v). By integrating d along v, we have a(0,v). 

We need more information for h, /, W 7 , w, a and g at the u = surface. We obtain /i from Equa- 
tion pTl) . / from Equation P2^) . 14^ from Equation P5]) . w from Equation P51) . a from Equation (|31|). 
and g from Equation (|45[) . This finishes the task of assignment of initial conditions. 

We choose ro = 10, e = 0.3, and u = 0, leaving the two parameters (A, <5) free, where A is the 
amplitude of the pulse of the scalar field and 5 is the phase of the complex scalar field. 

3 Gravitational collapse in f(R) gravity 

3.1 Models and the /(i?)-induced singularity 

In this paper, we study the following two models. 
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3.1.1 Model 1: Starobinsky model 

We first use the so-called Starobinsky model }18) : 

f(R) = R + XR + (55) 

where A, i?o, and n are free parameters. In this paper, we choose A = 2, Rq — 0.001, and n = 2 
(Figure [T]). We can see that there is a stable equilibrium around $o — 0.9932 and Vb ~ 0.001966. 
Therefore, we can choose these initial conditions for $o and Vq. Then the other initial parameters are 
determined by consistency condition. In this setup, the cosmological horizon is I — y'M^o/Vo ~ 55. 
Therefore, for convenience, we should choose the size of black holes less than I ~ 55. 

It should be noted that if R is less than Rq = 0.001, then the Brans-Dicke field $ = f'(R) should 
be less than zero (upper left of Figured]). Therefore, the region of validity for R is i?o < R < oo. 
On the other hand, in this region, the Brans-Dicke field $ varies from zero to one. In other words, 
$ cannot be greater than one. This particular behavior is due to the choice of f(R): f(R) ~ R in 
the large R limit. 
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Therefore, in the large R limit, we observe that 

f(R) ~ i?-Ai? + Ai?2»+i_L + ..., ( 56 ) 

/'(i?) ~ l-2n\Rl n+1 ^± TT + ..., (57) 

/"(i?) ~ 2n(2n+l)A<> +1 -^ + ..., (58) 

Km J7($) = lim ^ (*~^r _ 2F($)^ d$ (59) 

= lim [ R (2f(R)-Rf'(R))f"(R)dR (60) 

« J W^T- ( 61 ) 

Therefore, if n > 0, then the effective potential J7($) in the $ = 1 limit is finite, although the 
effective force term 2f - Rf = W - 2V is infinite. 

Now let us observe the behavior of potentials V(&) and [/($) (lower left and lower right of 
Figure [TJ. Here, we plot $ < 1. If the Brans-Dicke field $ becomes much smaller than the local 
minimum in terms of the effective potential E/($) (for example, if the field rolls up to $ ~ 0.5), then 
it will roll down and oscillate around the local minimum. However, it is potentially possible for $ to 
touch the point $ = 1; then the region should be identified with the curvature singularity R = oo. 
This is a kind of singularity, but not from the point of view of general relativity. Such kind of 
singularity highly depends on the choice of f(R). We call this kind of singularity the f(R) -induced 
singularity. Such a singularity will happen for various f(R) models, depending on the particular 
case; perhaps, a suitable choice of f(R) may remove such singularities. 

3.1.2 Model 2: R+ (l/2)ci? 2 

In order to avoid the /(i?)-induced singularity and to study the general feature of the f(R) gravity, 
we choose a simpler model that has no /(i?)-induced singularity. Here, we study the simplest f(R) 
model: 

f(R) = R+^cR 2 . (62) 
Here, c is an arbitrary constant. Then 

/'(V>) = l + c^ = $ (63) 

and 

$ - 1 

i/> = . (64) 
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Therefore, 



and 



V($) = -fW + ^W-iP (65) 

= Yc^' 1)2 (66) 



U(t>) = J (*^-2V)d* (67) 

= Yc^~ 1)2 ' (68) 



For stability, we choose c > 0. 



3.2 Model 1: Basic features 

3.2.1 The simplest case: de Sitter space 

As a triviality test, we investigate the A = case. Then we observe the de Sitter space without a 
black hole (upper left of Figure [2]). The cosmological horizon is at r = I ~ 55 and parallel to the 
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v 

Figure 3: $ of the A = 0.25 and 6 = case (Model 1). 

in-going null direction. Despite being a basic result, this shows a non-trivial consistency check for 
our simulations. 

3.2.2 Neutral limit 

We now observe a neutral black hole and a less charged black hole (upper right of Figure [2]). First, 
we observe A = 0.25 and 6 — 0. In this case, as there is no phase difference between the real and the 
imaginary part of the complex scalar field <f>, the matter combination becomes neutral. Therefore, 
we can see just a space-like singularity and an apparent horizon. 

Figure [3] shows the response of the Brans-Dicke field during gravitational collapse. Basically, 
the Brans-Dicke field first moves towards the lower value instead of the asymptotic value, and then 
returns to a certain limit. The similar situations have already been studied by the authors in |15) . 
In that paper, it has been observed that the collapse of a neutral matter field induces a decrease 
of the Brans-Dicke field for oj > —3/2; and as oj approaches the —3/2 limit, the Brans-Dicke field 
becomes more and more sensitive. The similar behavior can be checked for f(R) gravity as well, 
since f(R) gravity is the oj = limit of the Brans-Dicke theory. Of course, there is one major 
difference between the present and the previous work. In the previous paper, there is no potential 
and hence the field value returns to the asymptotic value, while in the present paper, the field value 
is governed by a potential and hence whole processes is more complex. 

3.2.3 Near extreme limit: /(i?)-induced singularity 

Here, we observe A — 0.25 and 6 — O.hit (lower right of Figure [2]). In this limit, we can push as 
much charge as possible |14j . In the Einstein gravity, because of electromagnetic repulsion, it is 
impossible to observe a naked singularity or an exact extreme black hole [H] [13] [H] ■ 



14 



However, one interesting point is that the /(i?)-induced singularity is located outside the appar- 
ent horizon. Gravitational collapse will perturb the Brans-Dicke field value to become smaller than 
the local minimum of the effective potential !/(<£>). Then, the perturbed Brans-Dicke field will roll 
down and rolls up to approach $ — > 1. Figure 0] shows such a behavior around the /(i?)-induced 
singularity. Figure 0] shows the behavior of log |1 — <I>| as a function of v. As v increases, |1 — <f>| 
first increases (i.e., $ decreases), and then |1 — <I>| decreases (i.e., $ approaches one) and eventually 
the simulation breaks down. This explains the behavior of the /(-R)-induced singularity, since the 
$ — > 1 limit is the R — > oo limit for the model. This singularity eventually induces the Cauchy 
horizons along the in-going and the out-going directions. 

Figure [5] explains what is going on around the singularity. The upper diagram shows three terms 
for the Brans-Dicke field equation, which is 

2 

r$ =<£=-(fZ + gW) - ^- (8^T C + [/'($)) . (69) 

The first part of the right hand side (fZ + gW) is negative while the second part of the right hand 
side (ret 2 1 12) (8ttT c + U') changes its sign. The latter increases as v increases. This becomes an 
effective resistance for $ that prohibits it to touch $ = 1 (i.e., ~ [/'(<&)). However, $. ul , will be 
solely determined by the sum of the first and the second term, and if the conditions are finely tuned, 
the first term and the second term can be exactly canceled around $ ~ 1. Then 

< 1 and ~ < 1 70 

will hold. (In this calculation, the former is of the order of 10~ 3 and the latter is that of 10~ 2 .) 
Then $ becomes a slowly varying function; as v increases up to order of one, the Brans-Dicke field 
$ varies linearly (lower diagram of Figure[5]) as a function of v. In this case, the increased value A$ 
relative to the scale of $ (« 1) is of the order of 10 -2 , and hence all the arguments are consistent. 

The physical importance is not easy to interpret. Since the /(i?)-induced singularity can lie 
outside a black hole, we can infer that it may have some implications towards astrophysical obser- 
vations, if we assume f{R) gravity to be the correct theory for our universe. On the other hand, 
if we believe cosmic censorship, then such f(R) will not be allowed. This crucially depends on the 
choice of the f(R) function. Moreover, even though a model allows /(i?)-induced singularity, if we 
suitably choose the model parameters, then the specific f(R) model will not allow the existence of 
such naked /(i?)-induced singularity for astrophysical scales. Therefore, further careful studies are 
needed. 
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Figure 4: log |1 — $| as a function of v along u — 14, 15, 16, 17, 18 for the A = 0.25 and 5 = 0.5tt 
case (Model 1). 

3.2.4 Intermediate limit: Oscillating horizons 

Here, we observe A = 0.25 and 5 = 0.37T (upper left of Figured]). Now, we can see the /(i?)-induced 
singularity and the Cauchy horizons, and they are inside the outer apparent horizon. There are 
some oscillatory apparent horizons inside the outer apparent horizon. This is due to the Brans- 
Dicke field and the potential. Since the Brans-Dickc field oscillates, the energy-momentum tensor 
components oscillate, and hence the sign of r„ can oscillate. 

An interesting question is whether such oscillating horizons do converge to form an inner Cauchy 
horizon in the v «Max limitQ as in the case of the Einstein gravity or not. Here, t>Max is not infinite, 
since we consider the black hole in a de Sitter space. Figure [5] shows that r „ will oscillate, but will 
never converge to zero. Therefore, in the v — ¥ VMax limit, there will be a Cauchy horizon, but not 
an inner apparent horizon inside a charged black hole. 

Finally, we summarize various features of gravitational collapse of the f(R) charged black holes 
in Figure [7J We list the following basic and new features of dynamical charged black holes in f(R) 
gravity: 

• Since f(R) gravity is the w = limit of the Brans-Dicke theory (with a potential), the 

"Max i s the largest value of the advanced time v inside a charged black hole. For an asymptotically de Sitter 
space, i>Max should be a finite value. For an asymptotically flat space, t>Max can be infinite. 
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Figure 5: Some terms (ra 2 /12)(8irT c + U'), fZ + gW, and $ im) (upper) and $ (lower) as functions 
of v along u = 15 for the A = 0.25 and S — 0.5ir case (Model 1). We focus on the behavior near 
the singularity. 
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Figure 6: |r>| along u = 10, 10.5, 11 of the A = 0.25 and 5 = 0.3tt case (Model 1). Since |r>| is 
plotted on the log scale, one can see that r v does not tend to zero in v — > «Max limit. 
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gravitational collapse will thus push the Brans-Dicke field value to a lower one than the 
asymptotic value, and this causes the perturbation of the Brans-Dicke field. 

• Perturbations of the Brans-Dicke field can cause the /(i?)-induced singularity, which in turn 
will make in-going and out-going Cauchy horizons. 

• There is a space-like singularity and a separation between the space-like singularity and the 
outer apparent horizon in v — ¥ i>Max limit, as in the case of the Einstein gravity; however, the 
v —> «Max limit is not an inner apparent horizon rather an in-going null Cauchy horizon. 

3.3 Model 2: Black holes without /(i?)-induced singularity 

Due to the /(-R)-induced singularity, it is difficult to see the inner structure of a black hole. Also, it 
is important to check whether we can build a model without /(i?)- induced singularity or not. We 
find that the R + (l/2)ci? 2 model is a good toy model to observe such a behavior. In this paper, 
we fix c = 10~ 3 . 

3.3.1 Neutral limit 

First, we want to reproduce a neutral black hole from the Model 2 (top of Figure [8]). In this case, 
we can again see the space-like singularity and the outer apparent horizon. 

3.3.2 Near extreme limit 

If we choose S — 0.5ir, we can get a near extreme black hole for the Model 2 (middle of Figure |8]). 
Here, we see a drastic difference as compared to Model 1. We can observe a space-like singularity 
and an outer apparent horizon. In addition, there is a clear separation between the singularity and 
the outer horizon inside the black hole. This is quite similar to that of a near extreme charged 
black hole in the Einstein gravity [12] [13] [14] . One interesting difference is there are oscillatory r >v 
horizons near the singularity, due to the term coming from the Brans-Dicke field. In pure Einstein 
gravity with charged scalar fields, since there is no source to violate the null energy condition, we 
do not get any such behavior. 

3.3.3 Intermediate limit 

Finally, we observe the intermediate charged black hole case (bottom of Figure . Again, there 
is a space-like singularity and an outer apparent horizon. In this case, the oscillatory behavior 
of r >v horizons is much clear. One can easily compare the size of the outer apparent horizon for 
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Figure 8: Results of Model 2: the neutral limit (top), the extreme limit (middle), and the interme- 
diate limit (bottom). 



20 



each result. Due to the charge repulsion, the size of black holes become larger and larger, as one 
decreases the charge |14j . 

Therefore, as conclusion for the study of R+ (l/2)ci? 2 model, we can notice that the existence 
of the /(i?)-induced singularity is not a necessary property of the f(R) gravity. 

3.3.4 Question: When does the /(i?)-induced singularity appear? 

At this stage, the natural question that one might address is the following: when does the f(R)- 
induced singularity appear? In other words, what are the conditions that induce such a singularity? 
To give a rigorous answer for the question, we need to check various kinds of f(R) models and 
hence it goes beyond the scope of the present paper. However, we can comment on some of the 
basic intuitions. 

Let us assume f(R) = R + const. + 0(R^ 1 ) + ... for the R — > oo limit. This implies that the 
model approaches the Einstein theory in the R — ¥ oo limit and hence the Brans-Dicke field should 
be a constant in the limit: $ — > $o- Note that, there is a relation between R and <!> by V'(&) = R. 
Therefore, in the R ->■ oo limit, V'($ ) = R -> °o; at the same time, E/'($ ) = 2f(R) - Rf'(R) ~ 
R -> oo. Also, we can check that J7($ ) = f°°(2f(R) - Rf'(R))f"(R)dR oc ©(iT 1 ). Therefore, 
J7($o) is finite and the effective potential £/($) has a cusp at <I>o. Then, the observation of the 
/(i?)-induced singularity is possible. 

The existence of /(i?)-induced singularity is related to the existence of a cusp in the effective 
potential U(<f>) due to a certain f(R) model. We can conclude that, the existence of f(R) -induced 
singularity is possible if f(R) ~ R for the large R limit, although it is not a general behavior of the 
f(R) gravity. 

3.4 Mass inflation in f(R) gravity 

We can define the null geodesies (here, our convention is [u, v, 9, ip\) as 

I" = ^(0,1,0,0), (71) 
a 

n» = ^-(1,0,0,0), (72) 
a 

so that = n^n^ = and l^n^ = — 1. Then, any time-like geodesic t^ can be decomposed as 

e = al» + bn^ (73) 

with the constraint ah = 1/2. Hence, any time-like observer will see the local energy density to be 

T^tn v = \ (b 2 T uu + a 2 T vv ) + \t uv . (74) 
or or 
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Therefore, if for some case T uu /a 2 , T vv /a 2 , or T uv /a 2 increases exponentially, for general choices 
of a and b, we can sure that it is the signature for the evidence of mass inflation. 

Figure |9] shows T uu /a 2 (upper), T vv /a 2 (middle), and T uv /a 2 (lower) for the A = 0.25 and 
6 = 0.37T case of the Model 1. The former increases exponentially, while the latter does not show 
any such behavior. Therefore, it is quite reasonable to conclude from this study that f(R) gravity 
has mass inflation. 

However, for the same case, we can check the behavior of the Ricci scalar R (Figure ITU]) and we 
can easily notice that it is finite. This is due to the potential of the Brans-Dicke field. If the Brans- 
Dicke field can be confined around the local minimum, then the Ricci scalar R should be of the 
order of ~ 0.004. This is not entirely trivial; as we saw in the /(i?)-induced singularity case, 

the Ricci scalar can diverge via complex interactions between the matter and gravity. However, we 
see that mass inflation does not induce the Ricci scalar to diverge. If we compare the result with 
the Einstein case (Figure 16 in [13]), we see that it is quite a distinct property from that of the 
Einstein gravity. 

The Model 2 gives similar results. Figure fTTl shows T uu /a 2 (upper), \T vv \/a 2 (middle), and 
\T U v\/ol 2 (lower) on a log scale. The tendency to linearly increase shows the existence of mass 
inflation clearly. 

However, the Ricci scalar in Figure [JJ] is not affected by mass inflation. It oscillates around 
an equilibrium point. One important note is that, although the Ricci scalar oscillates around an 
equilibrium position along a fixed u, there is a tendency that it increases as u increases. Figure [T51 
shows some detailed results for \R\ along some fixed v. As u increases, the amplitude of oscillations 
of R tends to increase, but it remains sufficiently small unless it approaches the central singularity. 
The scale of Ricci scalar is determined by c, as 

R= V'(<f>) = -($- 1). (75) 
c 

This implies that as the higher curvature correction term becomes more and more dominant (as c 
increases), the Ricci scalar R will slowly increase. 

To summarize the findings for the Model 1 and Model 2, we can conclude the following three 
things: 

• There is mass inflation, such that an in-falling observer will measure exponential divergence of 
some of energy-momentum tensor components. Then, v — > t>Max limit may be a null curvature 
singularity, similar to the Einstein case. 

• The Ricci scalar is bounded by a certain value along the out-going null direction. 
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Figure 9: T uu /a 2 (upper), T vv /a 2 (middle), and T uv /a 2 (lower) of the A — 0.25 and S — 0.3n case 
(Model 1). T uu /a 2 is plotted on the log scale; this shows an exponential increase in a component 
of the energy-momentum tensor as ~ exp^w, where k is a constant. 
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Figure 10: Ricci scalar R for the A = 0.25 and 5 — 0.3ir case (Model 1). 

• The Ricci scalar increases along the in-going null direction; but, the in-going observer will 
quickly reach the space-like central singularity. 

4 Discussion 

In this paper, we studied the gravitational collapse and dynamical behavior of a charged f(R) black 
hole. To satisfy the cosmological stability criterion, we chose some special models of f(R): one of 
them is the Starobinsky model, the other being the simplest R + (l/2)aR 2 model. 

The higher curvature corrections will change the basic features of the dynamical charged black 
holes from that of the original Einstein gravity. Since f(R) gravity is the lj = limit of the Brans- 
Dicke theory (with a potential), the gravitational collapse will push the Brans-Dicke held to a lower 
value than the asymptotic one, and this causes the perturbation of the Brans-Dicke held. This was 
already known from the authors' previous work |15j : however, the perturbations of the Brans-Dicke 
field can cause a /(i?)-induced singularity, and this is indeed a new observation for the f(R) gravity. 
This /(i?)-induced singularity will in turn lead to in-going and out-going Cauchy horizons. There is 
a space-like singularity and a separation between the space-like singularity and the outer apparent 
horizon in the v — » t>Max limit, as in the case of Einstein gravity; however, the v — > VMax limit may 
not be an inner apparent horizon, and will be just an in-going null Cauchy horizon. 



24 



10 20 30 40 



50 



00 



v 





Figure 11: T uu /a 2 (upper), \T vv \/a 2 (middle), and \T uv \/a 2 (lower) of the A = 0.25 and 5 = 0.3tt 
case (Model 2). All figures are plotted on the log scale; this shows an exponential increase in a 
component of the energy-momentum tensor as ~ exp^u, where k is a constant. 
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Figure 12: Ricci scalar R for the A = 0.25 and 8 = 0.37T case (Model 2). 




Figure 13: Ricci scalar \R\ for the A = 0.25 and 8 — 0.3ir case (Model 2). The \R\ axis is the log 
scale. In this case, as u increases, there is a tendency of increase in the Ricci scalar. 
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The main motivation of this paper is to study mass inflation in f(R) gravity. We observed some 
interesting features. First, there is mass inflation and an in-falling observer will measure exponential 
divergence in some energy-momentum tensor components similar to the Einstein case. However, 
interestingly, the Ricci scalar is bounded by a certain value along the out-going null direction. 

An out-going observer will never see a geodesically incomplete central singularity; rather, the 
observer will just see a curvature singularity caused by mass inflation with non-zero area. What is 
indeed a new feature is that the Ricci scalar is bounded for the observer. This is due to the higher 
curvature correction terms of the f(R) gravity, and also is due to the stable local minimum of 
the effective potential £/($). Thus, we can conclude that cosmologically 'stable' 'higher curvature 
corrections' can indeed hold the Ricci scalar to be of a finite value, even in the presence of mass 
inflation. 

Let us emphasize our results again. The action of the present paper is 



S = I dx i y/—g 



+ ^matter + f2(Rp Upa ) + h{Rpv 



(76) 



We only use the first and the second leading terms, since we assume that contributions of /b, /3, ... 
are sufficiently smaller than the leading terms. Then, 

1. The f(R) sector is finite even in the presence of mass inflation. 

2. The £ ma tter sector will increase larger than the f(R) sector. This is easy to check: if there 
is a kinetic term of an arbitrary scalar field ip, then it will be proportional to 1/a 2 , where 
a exponentially decreases in the presence of mass inflation. This shows the instability of the 
matter sector of the Lagrangian. 

3. Some of the energy-momentum tensor components are related to some of the curvature compo- 
nents by the Einstein equation. Therefore, some terms of higher order curvature corrections 
/2,/3,--- can be larger and larger in the presence of mass inflation. This may violate the 
original assumption of the action in Equation ([TB")) . 

4. We can conjecture that if we equally consider all the higher order curvature corrections 
fi, fs, ... and if we impose the cosmological stability^ for the entire action, then all the curva- 
ture components can be finite even in the presence of mass inflation. 



2 Here, the meaning of the cosmological stability should be well-defined. If the general action Equation ( 17"t>l > is 
equivalent to the Einstein gravity with a number of scalar, vector, or tensor fields, then we can impose the stability 
to the scalar, vector, or tensor fields; they should be asymptotically stable up to field perturbations. However, in the 
present paper, we cannot justify whether we can impose the concept for the most general case. But, at least, we can 
say that this is a reasonable generalization. For this point, we thank to the anonymous referee. 
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5. Then, the observed energy density of an in- falling observer near the Cauchy horizon can be 
finite even in a charged black hole. This can be a resolution of the mass inflation curvature 
singularity. 

Note that 1, 2, and 3 are confirmed by the present paper, and 4 and 5 arc our conjecture. 

In conclusion, we have a promising, yet brave, conjecture, which is as follows: if we extend 
the action to include not only a function of the Ricci scalar, but also functions of more complex 
combinations of curvature, and if we require cosmological stability for the entire action, then all the 
curvature components and the observed local energy densities should be finite around the Cauchy 
horizons of charged black holes. Then, we can trust the original action, as all the components of 
the action may be bounded by a certain value, and we can obtain a self-consistent gravitational 
dynamics in the presence of inner Cauchy horizons. This will shed some light on the study of regular 
black holes. Classical regular black holes have inner horizons in general, and they cannot be made 
free from the mass inflation instability. However, if higher order curvature corrections can make all 
the curvature components to be bounded, then it will be possible to understand the consistency 
condition for the action of a certain model, and it will for sure help us to understand more realistic 
regular black hole models. 
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Figure 14: Consistency test for \RW - R^ | around u = 10, 10.5, 11, for the A = 0.25 and 

5 — 0.37T of the Starobinsky model. 

Appendix. Consistency and convergence tests 

In this appendix, we report on the convergence and consistency tests for our simulations. We used 
the A = 0.25 and 6 = 0.3-7T case of the Starobinsky model, since this case includes all the interesting 
features of the present paper. 

For consistency, we can check various relations, but the most important non-trivial test for the 
f(R) gravity will be the Ricci scalar. The definition of the Ricci scalar is found in Equation (l4"9"l) . 
while we calculated it by using Equations (|50D and (l5"Tj) . We call the former i?W and the latter R^ 2 \ 
and checked \R,W - R^\/\R^\ around u = 10,10.5,11, where mass inflation occurs. Figure 1141 
shows that the differences are less than 10 -8 % and hence it is sufficiently small. 

For convergence, we compared finer simulations: 1 x 1, 2 x 2, and 4x4 times finer for some 
important slices (during mass inflation and near the f(R) induced singularity). In Figure [T51 we see 
that the difference between the lxl and 2x2 times finer cases is 4 times the difference between the 
2x2 and 4x4 times finer cases, and thus our simulation converges to second order. The numerical 
error is < 10 _3 %. 
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Figure 15: Convergence test: |r (lxl ) —r^x2) \/ r (2x2) and 4|r( 2X 2) -7"(4x4) I A(4x4) do coincide around 
u = 10, 12.5, 15, for the A — 0.25 and S = 0.3ir of the Starobinsky model. This shows the second 
order convergence. 
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